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Abstract. By combining methods of kinetic and density functional theory, we 
present a description of molecular fluids which accounts for their microscopic structure 
and thermodynamic properties as well as for the hydrodynamic behavior. We focus 
on the evolution of the one particle phase space distribution, rather than on the 
evolution of the average particle density, which features in dynamic density functional 
theory. The resulting equation can be studied in two different physical limits: diffusive 
dynamics, typical of colloidal fluids without hydrodynamic interaction, where particles 
are subject to overdamped motion resulting from the coupling with a solvent at rest, 
and inertial dynamics, typical of molecular fluids. Finally, we propose an algorithm 
to solve numerically and efficiently the resulting kinetic equation by employing a 
discretization procedure analogous to the one used in the Lattice Boltzmann method. 



1. Introduction 

Understanding transport phenomena and their interplay with structural properties of 
liquids in strongly inhomogeneous systems is of capital importance in many fields 
such as microfiuidics, science of materials, physiology and has several technological 
applications, ranging from pharmaceuticals to ternary oil recovery [1]. Much progress 
has been achieved in the last thirty years concerning the equilibrium properties of 
liquids in the presence of confining geometries or structured substrates. One of the 
major tools to explore inhomogeneous classical fluids is the Density Functional Theory 
(DFT), yielding thermodynamic (free energy, phase coexistence, osmotic pressure) and 
structural properties (density and composition profiles) within the same framework 
[2, 3]. Such a method presents many advantages: it is highly versatile, variational, 
lends itself to approximations guided by physical insight and is universally applicable 
to the study the properties of all real fluids under a wide variety of thermodynamic and 
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geometric conditions. In view of the great success of the DFT, a Dynamic DFT (DDFT) 
has been proposed to describe the dynamics of the density profile n(r, t) possibly due 
to a time dependent external potential or to a sudden change of the control parameters 
[4, 5]. 

The typical DDFT equation is a continuity equation for the number density: 

^^ + V.J(M) = (1) 

supplemented by a prescription relating the particle flux, J, to the density itself obtained 
by using the specific properties of the dynamics of the system. In the case of a colloid- 
solvent system without hydrodynamic interaction, one ideally traces out the solvent 
microscopic degrees of freedom and mimics their influence on the solute particles via 
Brownian dynamics induced by a heat bath stochastic term. Under this assumption 
one can derive a simple relation between J and n. On the other hand, in the case 
of molecular fluids, the absence of the fast degrees of freedom of the solvent seems to 
preclude such a possibility. It is necessary to recognize that the minimal set of collective 
variables which are relevant in the dynamical description of a fluid corresponds to the 
locally conserved flelds, namely number of particles, momentum and energy density. One 
can represent their evolution at a mesoscopic level by means of a set of hydrodynamic 
equations, which require the speciflcation of the currents associated with momentum 
and the energy densities in terms of the hydrodynamic variables via the so called 
phenomenological constitutive relations. However, if one requires a microscopic level 
of description similar to the one achieved by DDFT it is necessary to choose a kinetic 
approach and consider the time evolution of the onc-particlc phase space distribution, 
/(r, V, t). Doing so it is also possible to bridge hydrodynamics and microscopic structural 
theories of fluids. The method avoids the difficulty of computing separately the velocity 
moments of the distribution function because it considers directly the evolution of /. 
The derived equations can then be solved numerically by means of the Lattice Boltzmann 
method, a well-known method widely used in computational fiuid dynamics [6]. The 
resulting theory is able to predict structure and thermodynamics with an accuracy 
comparable to that of DDFT, but has the advantage of describing transport properties 
more realistically [7, 8]. In addition we shall show that from the present kinetic approach 
it is possible to derive both the DDFT for colloids and the transport equation describing 
molecular fluids. The two descriptions depart when one specifies the interaction of the 
fluid with the heat bath. In the first case, the heat bath is fixed, breaks the Gahlei 
invariance of the system and damps the momentum and energy transport. In the second 
case, the internal degrees of freedom act as an heat bath co-moving with the fluid and 
restoring the local equilibrium. 

The present paper is organized as follows: in sec. 2 we shall briefly review the DDFT 
for colloid-solvent solution by considering its derivation from an approach based on the 
evolution of the one particle phase space distribution. In sec. 3 we switch to the case 
of a molecular liquid and obtain a simplified transport equation which accommodates 
within the same framework thermodynamic and structural properties, inertial effects of 
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the fluid (sound waves, shear modes) and transport coefficients. In sec. 4 we give a brief 
account of the numerical method used to solve the transport equation. Finally in sec. 
5 we present the conclusions. 



2. Microscopic description 

Let us consider first a system of great practical interest, N particles, of mass m and 
positions and velocities Vj, immersed in a solvent. Since the solvent particles are 
in general much lighter and smaller than the host particles, the degrees of freedom 
of the solvent can be eliminated by a applying a suitable coarse graining procedure. 
The resulting total force acting on each colloidal particle is the sum of the interactions 
with other colloidal particles, of the external fields and of a stochastic force represented 
as a Langevin process, which combines a frictional force proportional to the velocity 
with respect to the solvent with a stochastic white noise forcing + ^j. Since 

both terms are originated from the solvent their amplitudes are related by the following 
fluctuation-dissipation relation («)) = 2Y°^''mkBT'°^''SijS''^S{t - s), where the 

superscript " soy indicates the origin of the friction and of the thermostatting medium. 
The particles interact through a pair potential U and are subject to a force field Fext- 



dvj 
dt 



dvi 



m—- — 



dt 



m-fVi + ^i{t) 



(2) 



It is straightforward to derive the associated hnear Kramers- Fokker- Planck [9] 
evolution equation for the 6N dimensional phase-space probability density distribution 
/jv({r, v}, i), where {r, v} indicates a 6N dimensional phase space point: 

d 



dt 



+ 
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Notice that we switched from a description in terms of trajectories of the particles to a 
probabilistic representation in terms of f^. This can be justified in two possible ways: 
the probabilistic description is a) the result of an ensemble averaging of the trajectories 
over a noise ensemble, or b) follows from an averaging over initial conditions. Case a) 
apphes to the damped stochastic dynamics (Y"''^ > 0) typical of a colloidal suspension, 
while case b) applies to the Hamiltonian dynamics characteristic of an atomic liquid. 
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solv 



0, when eq. (3) reduces to the Liouville equation. 



The information contained in /at is fully microscopic since it describes the 
microstate of the system. However, one recognizes that after an initial many-body 
regime, whose duration is of the order of the duration of a collision event, such a 
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representation becomes redundant and is possible to contract the description from 6N 
dimensions to only 6 dimensions, that is, from the phase space distribution of N particles 
to the one particle phase space distribution, /(r,v,t). In this stage, termed kinetic 
regime, the one-body distribution relaxes towards a local Maxwellian and its evolution 
can be represented by the following equation: 

(| + V • V + • |-)/(r, V, t) = Q{v, V, t) + B^'^'^ir, v, t). (4) 

at m ov 

This equation contains the streaming terms in the l.h.s., whereas in the r.h.s. Q is the 
interaction term 

Q(r,v,i) = ^Vv- J dr'j civ72(r, v, r', v', t)VrC/(|r - r'|) (5) 
and B^"'''" represents the heat bath term: 

S-'^(r, v,i) = 7^°'1^^ ^ + ^ • v]/(r, v,t) (6) 

To proceed further one makes the approximation that the two particle distribution 
function, /2 ,can be expressed as a product of one particle distribution functions times 
the positional pair correlation function g: 

/2(r, V, r', v', t) ^ /(r, v, t)f{r', V, t)gir, r') (7) 

so that eq. (4) becomes a closed non-linear equation for the one particle distribution. 
This is a crucial assumption since it decouples the evolution of /(r, v, t) from the 
evolution of the higher order multiparticle distribution functions. Many-particle 
correlations are however retained through the structural information contained in the 
positional pair correlation function function g{r,r'). As an approximation, we shall 
assume that g{r,r') is a non local function of the profile n{r,t), depends on time only 
through the density profile and has the same form as in a nonuniform equilibrium state 
whose density is n{r,t). 

For 7^°^'' > 0, one can simplify further equation (4) and derive an equation involving 
only the ensemble averaged particle density, by using the fact momentum and energy of 
the colloidal particles are not conserved, so that the currents become rapidly "slaved" 
to the density, i.e. the evolution is completely determined in terms of n(r, t) and its 
derivatives. That occurs because the velocities of the particles rapidly relax towards the 
equilibrium distribution, in a time of order of the inverse friction time r^"'''" — 1/^^°^'". 
This is the case when the non-dimensional ratio vtt^°^'" /a << 1 where a is the typical 
size of the particles and vt the thermal velocity, since the particles undergo a small 
displacement in a time ~ The only relevant evolution on timescales larger than 

^soiv Tj-egaxds the spatial distribution of the particles. This is the reason why the DDFT 
gives a sufficiently accurate description of colloidal systems. 

To be concrete, let us define the velocity moments of /(r, v, t) by multiplying eq. 
(4) by 1, mv, m(v — u)^/2, respectively and integrating w.r.t. v, obtaining the continuity 
equation: 

atn(r,t) + V-(n(r,t)u(r,t)) = (8) 
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the momentum balance equation 

mn(r, t)[dtu^{v, t) + ii,(r, t)diU^{v, t)] + ^.P-f ^r, t) - F,(r)n(r, t) - (r, t) = hf{v, t) 

(9) 

and the kinetic balance equation 

hBn{v,t)[dt + u,{v,t)d,]T{v,t) + P,f )(r,t)a,«,(r,t) + a,gf )(r,t) - C(2)(r,t) = h^^\v,t) 

(10) 

where the local fluid density is 

n(r,t) = y dYf{r,Y,t) (11) 
the local fluid velocity, u, is 

n(r,t)u(r,t) = / c^vv/(r,v,t) (12) 



and the local fluid kinetic energy density is: 

e''^{r, t) = %n{r, t)T{r, t) = I civ(v - n)'f{r, v, t) (13) 



We have introduced the kinetic component of the pressure tensor 



P,f)(r,i) =m j dv/(r,v,i)(v-u),(v-u),- (14) 



and of the heat flux vector: 



?f^(r,t) = f / civ/(r,v,t)(v-u)2(v-u), (15) 

Notice that Pj;^^ and q^^'' cannot in general be expressed in terms of the hydrodynamic 
moments introduced so far. We also have two terms stemming from the interaction 

/I s:'pnonidealt^/„ 
dvQ{v,v,t)v, = --n{r,t)V, ' >^ = -djH-\r,t) (16) 

C^'\r,t) = ^1 dvQ(r,v,t)(v-ur (17) 

The second equality in (16) relates C^-^^ to a functional derivative of the non ideal part 
of the free energy, whereas the third equality connects it to a spatial derivative of the 
potential part of the stress tensor. It is important to notice that C^'^'^ and C'-^^ vanish 
in a uniform system. Finally we specify the terms stemming the coupling with the heat 
bath 

.<"(M) ^ ™/.v(„- „).B-(r.v.O ^ -™.-„(M)«.(M) (18) 
6(2) (r, t) = ^ J dv(v - uyB"'^''{r, v, t) = 7^°'^[3A;Br^'''^ - m < >]n(r, t) (19) 
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The presence in eqs. (9-10) determine a fast equilibration process of the 

momentum current nu and of the local temperature T towards their stationary values. 
This fact is used to contract the description from the five fields appearing in eqs. (8- 
10) to the single field appearing in the DDFT equation. The rigorous mathematical 
procedure to derive this result employs the multiple time scale analysis [10, 11, 12], 
exploiting the time scale separation between the zero mode associated with the density 
fluctuations and the remaining modes, including the momentum and energy fluctuations 
which are fast relaxing due to the friction with the solvent. Physically speaking, the 
heat bath renders possible the contraction from the phase space to the configurational 
space because it rapidly "washes out" any possible deviation of the distribution function 
from the Gaussian. The larger the friction the faster the restoring process towards the 
Maxwellian. We shall not go through the details of such a derivation which can be found 
in refs. [10, 12] and merely quote the final result: 

dtn{r, t) = —y2 Vi[kBTVin{v, t) - F,(r)n(r, t) - Cf\v, t)]. (20) 

i 

It is possible to compute inertial corrections to such an equation by employing 
systematically the multiple time scale analysis, which is tantamount of an inverse friction 
expansion, but this procedure does not help in the case of molecular fiuids where the 
inverse friction parameter 1/^°'''" diverges, and only internal dissipation mechanisms 
are at work. A salient feature of molecular liquids is that they support hydrodynamical 
modes, as a result of local conservation laws of particle number, momentum and energy, 
which are eventually damped due to internal friction and heat transport mechanisms. 
Hence it is necessary to apply a strategy which preserves the translational invariancc of 
the system and does not select a particular reference frame where the solvent is at rest 
[13, 14]. 

3. Molecular fiuids 

The key difference between colloidal systems and molecular fluids is that in the former 
the equilibration is due to the heat bath externally imposed by the solvent whilst in 
the latter equilibration is realized via viscosity and heat conduction among the same 
fluid elements. In molecular fluids typical relaxation times of hydrodynamic modes 
diverge as A;~^ for excitations of wave vector A; — ?> 0, whereas the relaxation times of 
non- hydrodynamic modes arc shorter and flnitc in the same limit. Our treatment of 
eq. (4) must retain such a feature while accounting for the structural and thermodynamic 
properties with an accuracy comparable with DDFT [15]. For = 0, we recast eq.(4) 
as: 

dtf{v, V, t) + V V/(r, V, t) + ■ — /(r, v, t) = /C(r, v, t) + i3^"*(r, v, t) 

(21) 



Dynamic density functional theory versus Kinetic theory of simple fluids 



7 



where 
and 

fiocir,v,t)^n{r,t)[ , ]3/2 gxpf- ^/^ ~ V (23) 

■' ^ ' ' ^ ^ ' ^47rfcBT(r,t)^ 2kBT{r,t)J ^ ' 

The two terms in the r.h.s. of eq.(21) represent a suitable decomposition of Q(r,v,t). 
The idea is to treat the contribution of Q from hydrodynamic modes separately from 
the non-hydro dynamic kinetic modes. The term K, guarantees that the local transfer 
of momentum and energy is not altered. The part B'"* is chosen so that it does not 
contribute explicitly neither to the evolution of the hydrodynamic modes nor to the 
structure of the fluid, i.e. it does not appear in eqs. (9-10). For this reason it can be 
treated within a simple relaxation approximation, as shown below. We have deliberately 
employed the same symbol B used for the solvent heat-bath, to stress the fact that even 
in this case a thermalizing mechanism for the hydrodynamic modes exists and stems 
from the interplay between these modes and the non-hydrodynamic modes. The role of 
B'"* is to reproduce in an approximate and simple fashion the fast relaxation process 
of / towards the local equilibrium distribution fioc- The following BGK relaxation 
approximation [16] for the non-hydrodynamic components of Q meets the requirements 
of leaving unaltered the balance equations 

^-*(r, v,t) ^ -i/o[/(r, v,t) - /zoc(r, v,t)] = -z/o5/(r, v, t) (24) 

where z/q is a microscopic relaxation rate of the order of Enskog collision rate, loe = 
Ay/Trna'^VT, where a is a molecular diameter and vt = ^/k^Tjm. Since by construction 
both /(r,v, t) and fioc{r,v,t) correspond to the same values of n,u,T, one obtains 
that bf^ and associated with B*"*, vanish. In conclusion, the balance equations 
for the hydrodynamic moments generated by eq.(21) with the ansatz (24) meet all the 
requirements imposed by the translational invariance of the system and thus possess the 
correct long-wavelength hydrodynamic behavior. 

Before particularizing the description we wish to comment on a similarity between 
the DDFT and the present approach. It was found in DDFT that the free energy always 
decreases with time: 

dJ-'[n] f wiT-r ^-^[^ 



dt 

In the present case instead we have that the entropy, defined as 



S[f]^-kB J drdYf{r,v,t)lnf{r,Y,t) 



increases with time > 0, which is a statement of the iJ-theorem. 



dS[f] 
dt — "^"^^^ 

Using the special form of B*"* it is possible to obtain perturbatively an 
approximation for 6f and the transport coefficients at low density. In the spirit of 
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the Chapman-Enskog method [17] wc insert the trial distribution fioc{i^, v, t) in the Lh.s 
of eq. (21) and neglect the term /C. Using eqs.(8) - (10) to eliminate the time derivatives 
of the fields we arrive at the following explicit representation of 6f 
xf/ .N 1 /ioc(r, v,t) r/m (v - u)2 5\ 

/ (v — u)^ \ "I 
+ m\^{vi- Ui){vk - Uk) Sikj diUk{r, i) J . (25) 

which is valid in the low density hmit and to first order in the gradients of u and T. By 
substituting 5f in eqs. (14) ,(15) and with the help of the definitions 

pip = -n(^) (y,uj + Vjui - ^ V • u) 

and 

we compute the associated kinetic contributions to the viscosity 



and to the heat conductivity 

AW = l%Bn. (27) 

In order to obtain the "potential" contribution to the transport coefficients one 
needs an explicit representation of C^^^ and Ci'^\ which can be studied either in the 
random phase approximation (RPA) [18] consisting in combining eqs. (5) and (7) or in 
the Revised Enskog Theory (RET), which apphes to hard sphere fiuids and to piece wise 
constant potentials. 

Within the RPA one finds C(2)(r, t) ^ and 

Ci'\r,t) = n{r,t)Fr\r,t) (28) 
where we have introduced the self consistent molecular field: 

Fr'(r,t) ^- j dv'n{v\t)g{v,v')^U{\v-v'\). (29) 

Concerning the RPA, we wish to comment that it can treat fluids characterized by 
slowly varying potentials, hke ultra-soft potentials used to model the steric repulsion 
between polymers or long range attractive interactions. It can describe qualitatively 
their multiphase behavior , but does not give a suitable representation of their transport 
properties in the high density region, as it is evident from the vanishing of Ci'^\r,t). 
In fact, the RPA treatment of correlations while accounting for thermodynamical 
properties, docs not contribute to their hydrodynamical properties. 

Sharp repulsive interactions cannot be accounted for as a molecular field because 
they imply strong correlations between the relative positions and the relative velocity. 
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On the contrary, one can obtain a satisfactory representation of both structural 
and transport properties of hard sphere fluids. A very accurate treatment of how the 
hard sphere interaction contributes to the evolution of the phase space distribution is 
provided by the RET [19] which approximates the coUisional integral as: 



{9hs{r, r - cra)f{r, v[, t)f{r - aa, v'^, t) - ghs{r, r + aa)f{r, vi, t)f{r + aa, V2, t)} 



where ghs{i^,i^') is the hard sphere pair correlation function, a is the hard sphere 
diameter. The primes on velocities denote scattered values after the collision and 
V12 = (vi — V2), a is the unit vector directed from particle 1 to particle 2, and 



The difference from the Boltzmann equation is that the two distribution functions 
are evaluated at different positions. Such a difference in positions at collision allows 
the instantaneous transfer of momentum and energy between two particles. This 
coUisional transfer is much faster and effective than the translational transfer and 
becomes dominant as density increases. 

The salient feature which differentiates RPA and RET is the presence of velocity 
correlations in the latter approximation. In fact while it is assumed that the velocities 
of two particles are uncorrelated immediately prior to collision, but are correlated after 
they collide, because collision itself generates correlation. 

From eq.(30) we can obtain C^^^ and C*^^^ once /(r,v,t) is known. Interestingly 
one derives a simple relation for C*^^-* in terms of the hydrodynamic variables: 



where ' and P^^- ' are the coUisional components of the heat flux vector and of the 
pressure tensor, respectively. The two terms in the previous formula account for the 
transfer of thermal energy via coUisional transfer and the viscous dissipation due to the 
presence of velocity gradients, but such a connection is valid only for the hard sphere 
case, where there is no potential energy contribution. 

However, it is also possible, as shown in ref. [8], to obtain explicitly these quantities 
by approximating in eq. (30) /(r,v,t) by /;oc(r, v,t), that is, by neglecting 8f in the 
integrals. This approximation allows to express C\^\Y^t) and C^^)(r,t) by means of 
appropriate convolution integrals involving finite spatial differences of n(r, t), u(r, t), 
and T{v,t). It is such a dependence on the hydrodynamic fields which gives a non 
trivial result for the transport coefficients. In fact, one finds that in the limit of slowly 
varying inhomogeneities the present theory reduces to to the Longuet-Higgins Pople 




(30) 



V'l = Vi - (O- • Vi2)ct 




C7(2)(r,i) = -V,qf\v,t) - P^P{r,t)V,Uj{r,t) 



(31) 
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theory of transport coefficients [20]. This level of approximation yields the following 
coUisional contributions to the shear viscosity 

V^^^ = ^V^mvT(T^n''gU(T) (32) 
and to the heat conductivity 

A(C) = V'-^^mvT<j'n''gns{<j). (33) 
6 m 

To conclude, in the HS case we obtain a formula for the shear viscosity by taking 
into account only the distortion of the singlet distribution function from the uniform 
Maxwellian. For more general interactions the potential contribution to 77 dominates 
over the one body contribution in liquids at high densities and one has to consider 
the distortion of the doublet distribution function which is completely neglected in the 
present treatment [21]. 



4. Numerical method in a nutshell 



This section is included to illustrate how the microscopic description presented above 
can be applied in practical situations where the system under scrutiny is spatially 
inhomogeneous. In such a case, the transport equation (21) still represents a formidable 
numerical problem for the evolution of /(r, v,t) on a six- dimensional space. A possible 
method of solution is represented by focusing on the set of partial differential equations 
for the hydrodynamic moments (8-10), but these require the specification of complex 
boundary conditions and constitutive relations. 

An appealing alternative is to solve directly eq.(21) by using a discretization of 
velocities on a mesh [22, 23, 24, 25], a technique known as Lattice Boltzmann Equation 
(LBE), which deals directly with the phase space distribution function and provides a 
well-tested and robust method of solution [6]. 

In order to illustrate the scheme, we first project the phase space distribution 
function over an orthonormal basis spanned by the tensorial Hermite polynomials 

00 ^ 

/(r, V, t) = a;(v) 5^ (r, t)if« (v) (34) 

1=0 ^ 

where u;(v) = (27rv|.)~^/^e~"^^/^^T' and the Hermite coefficients 0a'(r, i), using the 
orthonormal relation 

I a;(v)4')(v)i^("'\v)rfv = {vtY^-^Si^S^^ (35) 

can be obtained by: 

0fnr,i) = / /(r,v,t)4')(v)dv (36) 

The exact infinite series representation of /(r,v, t) is approximated by a function 
/(r, V, t) obtained by retaining in eq.(34) only terms up to / = K, so that the complete 
and truncated distributions have the same coefficients up to Hermite order K. At 
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this stage, a Gauss-Hermite quadrature formula is employed in order to evaluate the 
expansion coefficients, 0a'*(r,t), recognizing that f{r,v,t)H^{v)/(jj{v) = p(r, v,t) is a 
polynomial in v with a degree < 2K, and using the nodes, Cp, and the weights Wp of a 
quadrature of order 2G > K, with p = 1, - ■ ■ ,Q. Using the result 

G 

dvuj{v)p{r, v,t) = Y^ Wpp{r, Cp, t) (37) 
one obtains 



where 



G 

r,t)= dvfir, V, t)Hl\v) = m^^Cp) (38) 

— Ill 

fp{r,t)^f{r,cp,t)-^ (39) 

UJ[Cp) 

Notice that fp{r, t) contains the same information as the continuous velocity distribution 
/(r, V, t), so that the continuous phase space distribution is replaced by a Q-dimensional 
array, /(r, v,t) — )■ fp{r,t). The Q nodes connect neighboring mesh points r on a 
lattice through the discrete velocities Cp, mirroring the hop of particles between mesh 
points, generally augmented by a null vector Cq accounting for particles at rest. The 
specific form of the lattice velocities and weights depends on the order of accuracy 
of the method and reflects the required Hermite order, as described in the following 
and thoroughly discussed in ref. [26] The LB algorithm exploits the Cartesian mesh 
to rearrange populations over spatial shifts corresponding to a first-order fully-explicit 
temporal update, 

where At is the time-step. 

Similarly, one expands the collisional and BGK contributions featuring in the right 
hand side of eq.21 and writes the following expressions: 

F-*(r, V, t) ^ S-*(r, t) ^ F"*(r, Cp, t)^ (41) 

and 

— II J 

/C(r, V, t) ^ /C^(r, t) = £(r, c^, t) (42) 
Using the specific form of the BGK term the LB algorithm reads 
fp{r + CpAt, t + At) = {l- —)fp{r, t) + —f;''{r, t) + [ICp{r, t) + Bp{r, t)]At (43) 

T T 

The resulting population dynamics of fp is able to virtually reproduce any target 
macroscopic evolutions to high accuracy. In absence of explicit collisional terms, 
and owing to the hyperbolic nature of the evolution equation, the Lattice Boltzmann 
approximates the Navier-Stokes equation in the nearly-incompressible limit. In addition, 
the LB method is completely flexible in terms of the mesh spacing Ax that can be 



Dynamic density functional theory versus Kinetic theory of simple fluids 



12 



tuned in order to resolve the details of the microscopic interactions, as for the spatial 
convolution in the coUisional terms. Therefore, the error introduced in the spatial 
discretization does not represent a critical issue. 

The specific form of the discrete velocities and weights depends on the order of 
accuracy of the method [26]. They are also designed to preserve mass, momentum 
and local isotropy by satisfying the rules Ylp'^p^pi — '^p'^p^pi^pj — '^T^ij-i 
'}2pWpCpiCpjCpk = and J2pWpCpiCpjCpkCpi = v^{6ik6ji + Sudjk), where Wp is a set of 
normalized weights and vt — i/ y/S is the mesh sound speed. 

The exact local equilibrium distribution function contains an infinite number of 
terms in its Hermite representation, but since we consider a truncation to second order 
Hermite polynomials we have 



= Wpn{r, t) 



Cpj'iij(r, ^ {cpiCpj 'yjn(5jj)'Uj(r, tjiij (r, t) 



(44) 



corresponding to a low-Mach [0[Ma^]) expansion of the local Maxwellian. On the 
other hand, by keeping only the momentum component of the coUisional contribution, 
C^^\r,t) , one writes the following discretized form: 



m 



Cp,C^^ (r, t) jcpjCpj - vlSij)ui{r, t)C^'^ (r, t) 

O ~l~ A 



(45) 



Similarly, the contribution arising from heat conduction, C^^^(r, f), could be included 
within the same framework. However, a thermal version of the LB scheme requires 
inclusion up to third or fourth order Hermite polynomials. 

Once the populations fp are known, they are used to compute hydrodynamic 
moments, entering both the equilibrium and in sampling the macroscopic evolution. 
The fiuid density, momentum current and local temperature read 

p 

n(r,t)u(r,t) = ^fp{r,t)Cp 

p 

^n{r,t)T{r,t)=J2fpi^M (46) 

p 

A Chapman-Enskog analysis shows that the relaxation time r is related to the kinetic 
component of the shear viscosity via 

ry(^) = nm4(T - — ) (47) 

i.e. the physical value, found in eq.(26), subtracted by a contribution of numerical 
origin. 

For hard spheres, the radial distribution function g{x^ r') appearing in the evaluation 
of Cj-^''(r, t) and C*^^''(r, t) associated with the RET approximation is obtained according 
to the Fischer and Methfessel construction [27]. At first, one defines a coarse-grained 
density n(r, t) via a uniform smearing over a sphere of radius (t/2, and the coarse-grained 
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packing fraction is r]{r,t) = 7ia^n{r,t)/6. The actual pair correlation, ghs{r,r + a), is 
replaced by its equilibrium value at the given smoothed density 

Qhsir, r + a) ~ [1 - 77(r + a/2)/2 + f (r + a/2)/4]/[l - r/(r + a/2)f. 

In order to evaluate the surface integrals, we choose a to be an even multiple of the 
lattice spacing and employ a 18-point quadrature over a spherical surface [28]. With this 
choice, the elements arising from g^g and the hydrodynamic moments are taken from 
6 on-lattice quadrature points while the elements arising from the remaining 12 off- 
lattice points are constructed via a linear interpolation from the surrounding on-lattice 
elements. 

The numerical method has been tested to determine the shear viscosity of a uniform 
system and an excellent agreement was found with theoretical prediction eq.(32). In 
addition we have considered the flow of an hard sphere system induced in a slit by the 
presence of a uniform field parallel to the walls. The observed velocity profiles had the 
expected parabolic behavior for wall distances sufficiently larger than the hard sphere 
diameter and moderately low densities, but displayed pronounced oscillations for larger 
densities and the average streaming velocity decreased as a consequence of mutual steric 
hindrance among particles. 

Regarding the stability of the method, LB is subjected to numerical instability 
whenever the flow velocity becomes larger than a certain threshold, function of the 
relaxation time, the wave number [26] and the stiffness of intermolecular forces. 
However, the stability range can always be widened by taking a smaller Ax. In 
fact, a generic upper bound for the variation of populations due to the forcing term 
is Sf/f ~ fp/wp ~ At/Cp < 1. 

The advantages of the LB method over the solution of the coupled hydrodynamic 
equations are evident: a) the convective term u(r, t) ■ V does not need to be computed 
explicitly; b) the solution of the Poisson equation for the pressure is not needed; c) 
the boundary conditions are handled in an easier way, for instance by imposing simple 
collision rules on the local velocity the distribution at the boundary, such as bounce- 
back in order to reproduce a no-slip condition; d) a relatively small number of velocity 
and spatial mesh points arc sufficient for many purposes; e) the two-body intcrparticle 
potential can be quite arbitrary and in the case of a slowly varying attractive tail it 
can be modeled via the RPA. Such an approximation, although not contributing to the 
shear viscosity, it allows to study phase coexistence. 

5. Conclusions 

Starting from a microscopic level we have obtained a governing equation for /(r, v,t) 
describing both equilibrium structural properties and transport properties. The essential 
difference between DDFT and the present method consists in the fact in the first the 
dynamical evolution of the fiuid system is treated only in configurational space and the 
fundamental variable is the one particle density, whereas in the latter the dynamics 
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is studied in phase space and the corresponding quantity is the phase space density 
distribution. The use of such a distribution enables us to consider the microscopic 
details of the molecular collisions and analyze transport properties on a quantitative 
basis. 

The DDFT equation can be obtained as a particular case, when the heat-bath 
term is externally imposed and has a viscous character forcing the velocity distribution 
to attain its global equilibrium form in a time much shorter than the time needed to 
reach the configurational equilibrium. As a result, DDFT does not describes sound 
waves in liquids, but only diffusive modes. In the case of molecular fluids instead, the 
time scales characterizing the equilibration times of the hydrodynamic flelds can be 
comparable and diverge as the wavelength of the modes tends to infinity, as a result of 
the translational invariance and absence of a fixed solvent. The kinetic equation (21) by 
construction reproduces this aspect and the mechanism restoring the global equilibrium 
in the fluid is purely due to the presence of gradients in the hydrodynamic flelds and 
not to the flelds themselves. Hydrodynamic versus non-hydrodynamic modes splitting 
proves a convenient route. As far as similarities are concerned, both methods use the 
same equilibrium structural information as their key ingredient and the same adiabatic 
treatment of two particle correlation necessary to decouple the one particle from the 
multiparticle dynamics. 

Other popular approaches within the LB framework have been put forward to treat 
inhomogcncous fluids, but they differ from ours in the derivation of the interaction terms. 
Most of these methods, such as the Swift-Osborne- Yeomans [29] and the effective pseudo 
potential methods by Shan-Chen [22], do not utilize a microscopic, but rather start from 
a coarse grained level and assume the validity of a local thermodynamics. In our opinion, 
we flnd rewarding to be able to derive in unifled way both DDFT and transport equation 
and to discriminate them on the basis of the thermostatting mechanism. 

The major limit of the present treatment is the assumption underlying the 
Boltzmann equation that particles are uncorrelated before collisions, resulting in a 
markovian dynamics. In the case of hard spheres, the RET approximation represents 
an improvement over the Boltzmann treatment by taking into account the non local 
character of the momentum and energy transfer, but without incorporating the time 
correlations which are responsible for memory effects. Even the RET approach is local 
in time with no memory of the past history. The deflciency of the theory becomes 
more acute when treating soft potentials, where the interaction is not instantaneous 
as in the hard-sphere case. For continuous interaction, it is not possible to consider 
the dynamics as composed of two separate stages, namely instantaneous collisions 
followed by free streaming trajectories, as in the case of hard spheres. The interaction 
have a flnitc duration and a second interaction can take place before the flrst one is 
completed. Some authors have modeled the different regions of the attractive potential 
tails by means of sequence of square- well constant potentials so that the RET formalism 
could be applied [30]. However, this approach inherits from the hard-sphere model 
the impulsive character of the interactions and therefore fails to reproduce memory 
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effects. Alternatively, it has been observed [31] that the motion of the particles can 
be regarded as the combination of continuous momentum changes, assimilated to a 
Brownian motion in an average attractive field, and quasi- instantaneous collisions with 
large momentum binary exchanges, described as hard core collisions. Following these 
lines it seems possible to extend the present approach to the treatment of soft potentials. 
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